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Demineralization of tooth mineral in 
the caries process was studied using a 
computer model that simulates a diffu- 
sion controlled dissolution process. The 
model consists of a two-compartment 
system. An acidic solution in the outer 
("plaque") compartment was assumed 
to be large in volume so that its com- 
position remained constant during the 
process. The solution in the inner ("le- 
sion") compartment was in equilibrium 
with the tooth mineral, but its composi- 
tion changed in response to diffusion of 
ions between the two solutions through 
an infinitely thin barrier. The permse- 
lectivity of the diffusion barrier to 
cations and anions can be modified as 
desired thus allowing the effects of 
membrane on the diffusion-dissolution 
process to be examined. Because the 
losses of calcium (Ca) and phosphate 
(P) from the "lesion" to the "plaque" 
generally does not occur at a molar ra- 
tio of 5/3, the Ca to P ratio of the dis- 



solving mineral, the composition of the 
"lesion" fluid can change significantly 
from the starting composition, and this 
in turn modifies the Ca and P fluxes. A 
steady state condition is eventually 
reached under which the ratio of flux 
of Ca to that of P becomes 5/3. The 
results of the simulation show that for 
a given "plaque" pH, the rate of dem- 
ineralization at steady state was the 
highest for cation and the lowest for 
anion permsclcctive membranes. These 
results were in good agreement with 
those from an experimental study under 
comparable conditions. 
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1. Introduction 

In the early stages of tooth decay, the affected 
tooth enamel usually is in the form of a subsurface 
lesion, Fig. 1. The outstanding features of a lesion 
are: (1) the "intact layer," a layer of relatively 
sound enamel at the tooth surface, (2) the "body of 
the lesion" consisting of partially demineralized 
enamel, and (3) the "advancing front" where active 
demineralzation occurs. The subsurface demineral- 
ization of tooth in the caries process requires diffu- 
sion of acid ions into, and of solubilized mineral 
ions out of the lesion. Even when highly simplified, 
it is a relatively complex physicochemical process. 



Results from recent studies have shown that during 
the demineralization process, the solution within 
the lesion is approximately saturated with respect to 
the tooth mineral at all times [1,2]. These results in- 
dicate that the rate of dissolution of tooth mineral 
in the lesion is faster than the rate of transport of 
the mineral ions out of the lesion. Thus, the process 
of subsurface demineralization is diffusion con- 
trolled, and the rate of lesion progression may be 
expected to be strongly dependent on factors that 
govern the diffusion process. 



593 



Volume 96, Number 5, September-October 1991 

Journal of Research of the National Institute of Standards and Technology 




LESION 



TOOTH 



(1) INTACT LAYER 

(2) BODY OF LESION 

(3) ADVANCING FRONT 

(4) SOUND ENAMEL 

(5) DENTAL PLAQUE 

(6) PELLICLE 




Figure 1. A schematic drawing of a typical dental caries lesion siiowing (1) an intact layer at tiie surface, (2) body of tiie lesion 
consisting of partially demineralized enamel, (3) tiie "advancing front" wliere active demineralization occurs, (4) sound enamel, (5) 
dental plaque covering tiie tootii, and (6) pellicle, adsorbed proteinaccous material on tiie tootii surface. 



Several mathematical models have been pro- 
posed to describe the diffusion process during 
lesion formation [3-7]. In the present study the 
caries process was also assumed to be controlled by 
the diffusion of ions, but it takes into consideration 
a most important phenomenon which has not been 
addressed in previous studies: the interaction be- 
tween the permselective diffusion of ions and the 
dissolution of tooth mineral as described below. 

The dissolution of tooth mineral at the advanc- 
ing front (Fig. 1) adds calcium (Ca) and phosphate 
(P) ions to the solution within the lesion at a molar 
ratio of 5/3, the Ca/P ratio of the solid. On the 
other hand, the relative rates of the loss of Ca and 
P from the lesion to the plaque generally is not at a 
ratio of 5/3. This is because diffusion of the various 
Ca and P containing ions (Table 1) is largely con- 
trolled by factors such as the electrochemical po- 
tential gradients and the permselectivity of the 
diffusion barrier, which may consist of the body of 
the lesion, the intact layer, and the pellicle, etc. [8]. 
The unequal rates of addition and removal of Ca or 



P by the dissolution and diffusion processes, re- 
spectively, would lead to a change in the composi- 
tion of the lesion fluid. This in turn modifies the 
driving forces for diffusion of all the ions. Thus, 
subsurface demineralization may be described as a 
coupled dissolution-diffusion process in which the 
composition of the saturated solution in the lesion 
changes in response to the diffusion of ions until a 
steady state is reached. 

In previous studies [9-11] a diffusion cell com- 
prising two compartments separated by an artificial 
membrane of known ion permselectivity was used 
as an experimental model to explore factors that 
may affect caries formation. One compartment 
(the "lesion") contained an excess of hydroxyap- 
atite (OHAp), Ca5(P04)30H, crystals, and its solu- 
tion was kept at or near saturation by stirring. An 
undersaturated acidic calcium phosphate solution 
flowed continuously through the other compart- 
ment (the "plaque"), thus providing the driving 
force for dissolution of the crystals as modified by 
the permselectivity of the membrane. As described 
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above, when the steady state is reached the compo- 
sition of the "lesion" solution will become constant 
and the ratio of the flux of Ca to P will be 5/3. The 
results from the bench-scale caries model studies 
demonstrated that the composition of the satu- 
rated solution in the "lesion" can undergo great 
changes prior to attainment of steady state [9,10]. 
Furthermore, the large changes in the concentra- 
tions of the various ions in the "lesion" fluid, as 
affected by the permselectivity of the diffusion bar- 
rier, significantly increased or decreased the rate of 
mineral loss at steady state [11]. In the present pa- 
per we describe a mathematical model that simu- 
lates this bench-scale caries model so that the 
dynamics of diffusion-dissolution interactions may 
be more fully understood. The model also makes it 
possible to more rapidly survey the effects of a 
number of different variables on demineralization. 

Table 1. Diffusing species and diffusion coefficients at 25 °C 



Species 


D^xlffcK 


i7s 


Species 


D' 


"xlO^cmVs 


H* 


9.3» 




OH- 




5.28" 


K+ 


1.96" 




ci- 




2.03" 


Ca^* 


0.79" 




H2PO4" 




0.88" 


CaHjPO; 


0.62° 




CaP04" 




0.29' 


CaOH+ 


1.03" 




HPOr 




0.44" 


H3PO4 


1.61'' 




PO4" 




0.61" 


CaHP04 


0.28= 











" Ds calculated from values given in [12]. 
"Ds calculated from values given in [13]. 
'D'l calculated from values given in [12] and [13]. 



2. Methods 
2.1 Overview 

In this computer simulation of the caries-forming 
process, the transport of all solution species (Table 
1) is tracked between a "plaque" compartment and 
a "lesion" compartment, separated by a barrier 
which blocks the free flow of fluids and is selec- 
tively permeable to all solution species. The 
"plaque" compartment is treated as a large reser- 
voir of invariant composition with a fixed degree of 
under-saturation with respect to the dissolving 
mineral, OHAp. This composition can be selected 
to simulate various oral environments. The 
"lesion" compartment is visualized as a suspension 
of OHAp crystals maintained at or near saturation 
by dissolution or precipitation which is assumed to 
be more rapid than the diffusion process. 

The computations proceed to convergence 
through a number of cycles, each consisting of two 



phases or loops. In the first phase ("diffusion 
loop"), diffusion across the barrier for a specified 
small time interval results in concentration changes 
in the "lesion" compartment. In the second phase 
("dissolution loop"), saturation of the "lesion" so- 
lution is re-established. When the system has con- 
verged to a steady state, the composition of the 
"lesion" compartment becomes invariant with time, 
and concomitantly the molar ratio of the total flux 
of calcium /(Ca)T to the total flux of phosphate 
/(P)t attains the value 5/3, the Ca/P of the dissolv- 
ing solid. 

The model approximates continuous diffusion 
across an infinitely thin barrier separating the "le- 
sion" from the "plaque." Under this restriction the 
diffusion potential (A£) and the fluxes (/s) of the 
individual solution species are approximated [14] by 
the following equations: 



A£ = -(RT/F)%,(tJzs)\n(a^/al) 



(1) 



and 



/s=-ACs[ln(as7fls) + Zs(F/RT)AE] (2) 



where 

T = Temperature (K) 

R =The gas constant (8.314x10^ 

J-mol-'-K"') 
F = Faraday constant (9.649x10" 

J-V-'-mol-') 

and for the sth species, with the compartment desig- 
nated by superscript, 

ts = transference number 

Zs = valence 

Ds = diffusion coefficient within the barrier 
(cmVs) 

Cs = (CV+Cs)/2 (mol/L); the average concentra- 
tion across the membrane 

fls = thermodynamic activity (mol/L). 

It may be observed that the flux, /$, of any 
charged species is influenced by both chemical and 
electrical (membrane) potentials as indicated by the 
two terms in Eq. (2). 

2.2 Characterization of the Barrier 

For the barrier to exhibit the permselectivity of 
tooth enamel to individual ions, the diffusion 
contants, D°, can be modified by arbitrary scaling 
factors. For example, to simulate an anion perm- 
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selective membrane, the diffusion constants of all 
cations may be reduced by a common scale factor, 
f, using the equation, 

D.=f\^^-Dl (3) 

where Dl is the diffusion coefficient of species s in 
water (Table 1). Lfi for the ions was calculated 
from the ionic mobility (|Xs) using the Nernst-Ein- 
stein relation \pi = {RTIF)\x.^ and for ion pairs, the 
Nernst method for the average diffusion coefficient 
[12]. Values of/ of 1.0 (nonselective membrane) or 
0.01 for either cations (anion permselective mem- 
brane) or anions (cation permselective mem- 
branes) were used in this study. The diffusion 
coefficients for H* and OH" were allowed to be 
invariant to account for the Grotthuss mechanism 
of H"^ and OH" association with the solvent (H2O) 
across the membrane [12]. 

2.2.1 Input and Initial conditions 



1. 



2. 



3. 



All constants in the preceding list from R (gas 
constant) to D° (diffusion coefficients) were 
stored as functions of temperature in the pro- 
gram. Also stored were dissociation and associ- 
ation constants of all partially dissociated 
species such as weak acids, ion pairs, etc. 
The temperature, T, was input for each run and 
was the same for both the "lesion" and 
"plaque" compartments. 
For the "plaque" compartment, the composi- 
tion of the solution was chosen to represent a 
fixed degree of undersaturation with respect to 
tooth mineral, GHAp, at the given tempera- 
ture. This solution retained these initial con- 
centrations. 

Table 2 lists the compositions of five 
"plaque-saliva" compositions that had been 
used in the previous experimental study and 



were used in the present simulation study. The 
solutions all contained 1 mmol/L HCl, 30 mmol/ 
L KCl, with various amounts GHAp dissolved 
such that the pH of the solution ranged from 
3.2 to 5.0 and pWPoHAp {LAP = \on activity 
product) ranged from 82.55 to 65.1. 

4. For the initial solution in the "lesion" compart- 
ment, the composition was chosen to be the 
"plaque" solution that was made saturated with 
respect to GHAp (Table 2). This solution can 
be represented as a point on the p[/(-)' 
isotherm, e.g., point A, in the phase diagram 
shown in Fig. 2. The generalized fourth compo- 
nent [[/(±)] in the system is defined in Ap- 
pendix A. 

5. Compositions of solutions in both compart- 
ments were initially adjusted to assure elec- 
troneutrality by an iteration procedure that 
finds the appropriate [H*] value (see Sec. 
2.2.2,1 below). 

2.2.2 Diffusion Loop of Calculation Objective: 
Compute state of the system in the "lesion" after an 
interval of mass flow. 

1. Concentrations and activities of all species Cs 
(except [H*]) listed in Table 1 were calculated 
by an iterative process [15,16] using known 
equilibrium constants and a suitable approxi- 
mation (Debye-Hiickel) for the activity coeffi- 
cients [15]. All constants were determined at 
the temperature of interest and are given in 
Table 3 for 25 °C. The computations were done 
just once for the "plaque" compartment be- 
cause this compartment is held at fixed concen- 
trations. 

2. Next the diffusion potential, A£, was calculated 
according to Eq. (1); A£ is a linear function of 
the terms consisting of logarithmic ratios of the 



Table 2. Composition and pMPoHAp of the "plaque" and initial lesion solutions 





pH 


Tot Ca 


Tot PO4 


K 


CI 


pMfoHAp" 






raraol/L 


mraoI/L 


raraoI/L 


mraol/L 




Plaque solutions 














A 


3.194 


0.741 


0.444 


30.0 


31.0 


82.55 


B 


3.503 


0.185 


0.111 


30.0 


31.0 


77.17 


C 


4.003 


0.630 


0.378 


30.0 


31.0 


72.48 


D 


4.503 


0.689 


0.413 


30.0 


31.0 


68.66 


E 


4.998 


0.709 


0.425 


30.0 


31.0 


65.10 


Initial lesion solution 














5.929 


0.741 


0.444 


30.0 


31.0 


58.55 



■pMPonAp values calculated at 25 °C, appropriate constants given in Table 3. 
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Figure 2. A schematic representation of the changes in the composition of the "lesion" 
solution, expressed in a -log[P]T vs -log[Ca]T phase diagram, during a computation 
cyele. See Appendix for the definition of U. 



Table 3. Constants at 25 °C 



Equilibrium constant (K) 



Reaction 



K 



Association constant (Ka) 



Ion pairs 



K 



^caHjPoj: Ca^'* + H2PO; jiCaHjPOJ 
ifcaHPoK Ca^'* + HPO3- ?± CaHPOj 
Kc,?or Ca^* + POt" ?± CaPOj 
Kc^oH*-- Ca^* + OH- ?± CaOH+ 



Solubility product (KSP) 



Ca-P04salt 



K 



KSPoHAf-- Ca5(P04)30H 

KSPocf. Ca8H2(P04)6-5H20 
KSPdctd: CaHP04-2H20 
KSPf,rcp: Ca3(P04)2 
KSPao-: Ca3(P04)2 



2.79x10-" 

3.71x10-"' 

2.51x10-' 

1.22x10-2' 

5.89X10-^ 



Ref. 



Kfi: H3P04fsH+ + H2POJ 


7.11x10-^ 


[17] 


Kp2: H2PO4 fs H+ + HPO5- 


6.31x10-" 


[18] 


Kfi: HPO5- ?± H+ + PO5- 


4.52x10-'^ 


[19] 


K„: H20f± H+ + OH- 


1.01x10-" 


[20] 


Debye (A) 


0.5115 


[21] 



Ref. 



8.45 


[22] 


265.0 


[22] 


2.90x10' 


[23] 


20.0 


[24] 



Ref. 



[25] 
[26] 
[22] 
[27] 
[28] 



activities in the two compartments, \n{a\/a^), 
and the transference numbers as shown in 
Eq. 1. 
3. Calculation of mass flows was conducted as fol- 
lows. 



a. The fluxes, /s, were computed using the 
modified Nernst-Planck flux equation 
[Eq. 2]; the approximations used in deriv- 
ing this equation as well as Eq. (1) for the 
diffusion potential require that all the 
quantities remain constant over a certain 
small time interval, A^ It was thus neces- 
sary to limit the magnitude of concentra- 
tion changes by adjusting the size of A^ 

b. The material flows are ACs=/sA^ and the 
resulting new concentrations are Cs' = 
Cs + ACs. These two quantities were used 
in the following constraints to assure that 
the changes in the state of the system were 
so small that the assumption in 3a holds. 

(i) |ACs|/Cs must lie between predeter- 
mined limits. 

(ii) For the new system composition, 
|A£„ew— A£oid| must not exceed a 
predetermined limit. 

If either (i) or (ii) was not satisfied, a new 
smaller At was chosen and ACs, etc., were 
recalculated. 

c. Using the new ion concentrations, Cs, new 
total concentrations of the various compo- 
nents resulting from the mass flow were 
calculated. The value of [Cajr for the 
calcium component is the sum of all 
calcium-containing species and similarly 
for P, K, CI and other components. The 
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system will now be at a point "D" in the 
phase diagram (Fig. 2), representing a con- 
dition of undersaturation. 
d. The total flux for each component was calcu- 
lated for later use to check steady state con- 
dition. In particular, for the Ca and P 
components: 

/(Ca)T=/(Ca^"')+/(CaH2P0t) + . . . 
/(P)t=/(H3P04) + . . .+/(CaH2P0j) + . . . 

where /(Ca)T and /(P)t are the sums of the flux for 
all species containing Ca and P, respectively. 

2.2.3 Dissolution Loop of Calculation The 

next stage of the calculation was to re-establish sat- 
uration in the "lesion" compartment. This feature 
is a crucial part of the model that has not been 
implemented successfully in previous literature re- 
ports. 

Objective: Drive the solution in the "lesion" 
from the undersaturated state (point D, Fig. 2) 
reached at (Sec. 2.2.2,3c) to a new equilibrium (sat- 
urated) condition. The position in the new compo- 
sition will then be at point B, lying on an isotherm 
pC/(-)"(Fig.2). 

Method: With the total phosphate concentra- 
tion, [P]t, and [H"^] as variables, a nonlinear least 
squares procedure [29,30] was used to satisfy the 
conditions of electroneutrality and saturation (Eqs. 
(4) and (5)), forcing congruent dissolution/precipi- 
tation, i.e., A[Ca]T = (5/3)A[P]T. 

SsZsCs = £ -> (4) 

pIAPoHAp-pKSPoHAp = 5^0 (5) 

where IAPqhap is the ionic activity product and 
^5PonAp is the solubility product for OHAp (Table 
3). 

The least squares procedure utilizes an iterative 
process to achieve convergence to a minimum 
where the sum of £2 + 5 2 approaches zero. The 
resulting composition is represented by a new 
point, B, on an isotherm pC/(-)" (Fig. 2). The 
isotherm pC/(-)" is different from the isotherm 
pC/(-)' in the previous cycle because of the 
change in the value of "electroneutrality unbal- 
ance" [8]. This is caused by the diffusion of non- 
consumed anions (e.g., CI") and cations (e.g., K"^) 
[10] that occurred during the time interval A^ The 
composition at point B was then compared to the 
composition at the commencement of the preced- 
ing mass flow stage (Point "A" in the diagram). If 
B differed from A, the steady state had not been 
reached and the system returned to step 2 in Sec. 



2.2.2 until convergence, i.e., the steady state was 
attained. 

2.2.4 Steady State At steady state the difference 
between points "A" and "B" in Fig. 2 becomes 
negligible. This implies that mass changes due to 
diffusion are exactly compensated by mass changes 
due to dissolution/precipitation reactions. Under 
these conditions /(Ca)T//(P)T = 5/3. The total 
fluxes were calculated at step 3d in Sec. 2.2.2 and 
the ratio was monitored throughout successive iter- 
ations. 



3. Results 

The composition of the "lesion" solution and the 
fluxes of the various components were found to 
change significantly from the starting point to the 
steady state. The rates of demineralization, ex- 
pressed in the units \ig OHAp/cm^ min, calculated 
from the steady state (Ca)T or the (P)t flux at vari- 
ous "plaque" pH values are given in Table 4. It can 
be seen that for a given type of membrane the rate 
of demineralization increased as the "plaque" pH 
decreased. For a given "plaque" pH the rate of 
demineralization is the highest with the cation 
permselective membrane and is the lowest with the 
anion permselective membrane. Although the ef- 
fects of membrane permselectivity on rate of dem- 
ineralization were similar in trend at the various 
"plaque" pH values, the effects became progres- 
sively greater as the "plaque" pH was decreased. 
For example, at pH 5.0 the relative rates of dem- 
ineralization for the three membranes were nonse- 
lective: cation: anion = 100: 241: 1.5, whereas at pH 
3.2 the relative rates were 100: 258: 0.064. For com- 
parison purposes the experimental values of Chow 
and Brown [11] are also given. It can be seen that 
the effects of pH and membrane permselectivity on 
the demineralization rate are similar between the 
computer simulation and the bench-scale experi- 
ments. The simulation and the experimental flux 
values were found to be strongly correlated with 
the correlation coefficients being equal to 0.995, 
0.996, and 0.999 for the nonselective, cationic 
permselective, and anionic permselective mem- 
brane systems, respectively. As shown in Table 4, 
the average ratios of the simulation and experi- 
mental flux values for the three groups are 29.0, 
43.0, and 1.3, repectively. The ratio, 29.0, for the 
nonselective membrane system may be taken as a 
proportionality constant needed to account for the 
inherent differences in the two models. The higher 
ratio, 43.0, for the cationic permselective mem- 
brane system and the lower ratio, 1.3, for the 
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anionic permselective system suggest that the simu- 
lation model magnified the effects of the mem- 
brane, probably due to the high degree of 
permselectivity assigned to the membranes (0.01). 
The effects of membrane permselectivity on the 
rate of demineralization shown in Table 4 are in 
fact a result of the significant changes in the com- 
position of the "lesion" solution that occurred from 

Table 4. Rate of demineralization obtained from the experi- 
mental bench-scale model and the computer simulation of caries 

Rate of demineralization of OHAp ((ig OHAp/cm^ min) 
"Plaque" pH 
Membrane 3.2 3.5 4.0 4.5 5.0 



Nonselective 














Simulation" 

Experimental'' 

Ratio'= 


265.3 
10.4 
25.5 


105.1 
3.08 
34.1 


28.11 
0.88 
31.9 


7.54 
0.28 
26.9 


2.66 
0.10 
26.6 


avg = 29.0 


Calionic 














Simulation 

Experimental 

Ratio 


683.7 
15.7 
43.5 


294.7 
8.08 
36.5 


79.64 
2.31 
34.5 


22.48 
0.62 
36.3 


6.42 
0.10 
64.2 


avg = 43.0 


Simulation 

Experimental 

Ratio 


0.17 
0.19 
0.89 


0.14 
0.14 
1.00 


0.09 
0.07 
1.29 


0.06 
0.03 
2.00 


0.04 
0.00 


avg = 1.30 



" Calculated from the flux of calcium or phosphate from the com- 
puter simulation. 

■■Experimental bench-scale data from [11]. 
° Ratio of the simulation and experimental flux values. 



the beginning of the process to the time when 
steady state was reached. For the purpose of illus- 
trating the dynamic nature of the diffusion-dissolu- 
tion interactions, results from the simulations with 
a "plaque" pH of 4.5 are described in greater de- 
tail below. 

Figure 3a shows the [Cajr, [P]t, [K]t, [C1]t, and 
the pH of the "lesion" solution as a function of 
time for the system with the nonselective mem- 
brane. The fluxes of these components and the 
membrane potential during the same time period 
are shown in Fig. 3b. It is noted that although a 
total of seven phosphate-containing species (Table 
1) were considered in the calculations, in the pH 
range used, H2PO4 is the predominant phosphate 
species and its diffusion accounts for the bulk of 
the P flux. Similarly, the diffusion of the free Ca^*, 
K^, and CI" ions accounts for nearly all of the Ca, 
K and CI fluxes, respectively. Consequently, the 
changes in the total concentrations and total fluxes 
shown in Figs. 3a and 3b can be adequately under- 
stood by examining the diffusion of the corre- 
sponding predominant ions. At the beginning of 
the process, driving forces existed for the inward 
diffusion of H* ions and outward diffusion of Ca^* 
and H2PO4 ions. The diffusion of these ions pro- 
duced a slight negative membrane potential 
(-0.07 mV), which induced a short burst of K* 
and Cr diffusion in opposite directions (Fig. 3b) 
even though no concentration gradient existed for 
either ion initially. The K* and Cl~ diffusion 
ceased as the concentrations of these ions in the 
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Figure 3a. The [Ca}r, [P}r, [K}r, and [Cl}r of the "lesion" solution as a func- 
tion of time for the system with the nonpermselective membrane and a 
"plaque" pH of 4.5. 
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Figure 3b. The fluxes of these components and the membrane potential over 
the same time period. 



"lesion" solution approached those prescribed by 
the membrane potential following the Nernst 
Equation [14]. In approaching the steady state, the 
[Cajr and [P]t increased slightly, from 0.741 and 
0.444 mmol/L to 0.855 and 0.513 mmol/L, 
respectively, and the pH decreased slightly, from 
5.929 to 5.855. The [K]t decreased slightly, from 
30.00 to 29.95 mmol/L, and the [C1]t increased 
slightly from 31.00 to 31.11 mmol/L. The (Ca)T and 
(P)t fluxes increased steadily to 1.324 and 0.795 
mmol/cm^ sec, respectively, when the steady state 
was reached. 

Figure 4a shows the [Cajj, [P]t, [K]t, [C1]t, and 
the pH of the "lesion" solution as a function of 
time for the system with the cation permselective 
membrane. The fluxes of these components and 
the membrane potential during the same time pe- 
riod are shown in Fig. 4b. It is seen from Fig. 4b 
that the cation permselective nature of the diffu- 
sion barrier allowed Ca^"^ ions to diffuse out of the 
"lesion" rapidly during the early period, i.e., the 
first 20 h of the process. Although the flux of H"^ 
was not followed in the calculations, the presence 
of a negative membrane potential at this time pe- 
riod would indicate that the inward diffusion of H* 
was in fact more rapid than the outward diffusion 
of Ca^"^. The negative membrane potential brought 



about a burst of outward diffusion of K"^, an ion 
which did not have a concentration gradient at the 
start. The membrane potential reached a minimum 
of -0.9 mV at approximately 50 h and started to 
increase steadily in approaching the steady state. 
The change in the sign of the membrane potential 
can be attributed to a decreased inward diffusion 
of H* because the "lesion" pH has been reduced 
from 5.93 to 5.45, diminishing the driving force for 
proton diffusion. The change in the sign of the 
membrane potential led to a corresponding change 
in the direction of K"^ diffusion. After the rear- 
rangement of dominant driving forces in the initial 
stage of the process (the first 200 h), the fluxes and 
concentrations of all ions approached the steady 
state values asymptotically. When the steady state 
was reached, the [Ca]T increased significantly from 
0.741 to 1.417 mmol/L, the [P]t increased substan- 
tially from 0.444 to 9.576 mmol/L, the [K]t in- 
creased from 30.00 to 34.61 mmol/L, and the [C1]t 
decreased from 31.00 to 27.70 mmol/L. It is inter- 
esting to note that at steady state the [Ca]T/[P]T 
ratio was 0.148, which is drastically reduced from 
the initial value of 1.67. The (Ca)T and (P)t fluxes 
at steady state were 3.719 and 2.164, nearly triple 
the corresponding fluxes in the system with the 
nonselective membrane. 
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Figure 4b. The fluxes of these components and the membrane potential over 
the same time period. 



Figure 5a shows the [Ca}r, [P]t, [K]t, [C1]t, and 
the pH of the "lesion" solution as a function of 
time for the system with the anion permselective 
membrane. The fluxes of these components and 
the membrane potential during the same time pe- 
riod are shown in Fig. 5b. It is seen from Fig. 5b 
that during the initial period, e.g., the first 5 h, the 
anion permselective nature of the membrane al- 
lowed the P ions to diffuse out readily. This and 
the inward diffusion of H"^ (as described in "Meth- 
ods," in the present model H"^ ion transport was 



not impeded by the anion permselective mem- 
brane) produced a negative membrane potential. 
The negative potential prompted an inward diffu- 
sion of Cr ions which gradually diminished as the 
steady state was approached. The changes in con- 
centrations of the various components were rela- 
tively small. At steady state the [Cajt increased 
slightly from 0.741 to 0.861 mmol/L, and the [P]t 
decreased from 0.444 to 0.425 mmol/L, resulting in 
an increase in [Ca]T/[P]T ratio from 1.67 to 2.03. 
The [K]t decreased from 30.00 to 29.90 mmol/L, 
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Figure 5a. The [Cajr, [PJr, [KJr, and [Cljr of the "lesion" solution as a func- 
tion of time for the system with the anion-permselective membrane and a 
"plaque" pH of 4.5. 
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Figure 5b. The fluxes of these components and the membrane potential over 
the same time period. 



and [C1]t increased from 31.00 to 31.17 mmol/L. 
The fluxes of (Ca)T and (P)t were 0.00963 and 
0.00665, respectively. These values are over 100 
times smaller that those in the system with the non- 
selective membrane. 

4. Discussion 

The mathematical caries model described in the 
present study is a simulation of the experimental 
bench-scale caries model reported previously 



[9-11]. The model was aimed at examining the in- 
teractions between the diffusion and dissolution 
processes that occur during the caries progression. 
Because in this model the solution species were al- 
lowed to diffuse across an infinitely thin barrier 
separating the "lesion" from the "plaque," the 
model is inherently limited in that it does not ad- 
dress certain characteristics of a true lesion such as 
the presence of concentration gradients within the 
lesion. The present model also does not describe 
the mechanism for the formation of the mineral 
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dense surface layer. On the other hand, the model 
clearly demonstrated the effects of membrane 
permselectivity on lesion composition and rate of 
demineralization. The results of the simulations 
parallel closely those obtained in the bench-scale 
caries model in several ways as described below. 

As described in Sec. 2, the compositions of the 
"plaque" solutions are such that they had the same 
HCl and KCl concentrations as those of the initial 
saturated solution in the "lesion" compartment but 
had lower amounts of OHAp dissolved (Table 3). 
Consequently, if one were to dissolve OHAp in any 
of the "plaque" solutions until saturation was 
reached, the solution would be the same as the ini- 
tial "lesion" solution. The fact that in many cases 
the steady state [Cajj and/or [P]t of the "lesion" 
solution became greater than the corresponding 
values in the initial "lesion" solution demonstrates 
the importance of the membrane effects on OHAp 
dissolution. 

The [Ca]T/[P]T ratio of the "lesion" solution at 
steady state was significantly lower than the initial 
ratio of 5/3 for the system with the cation-selective 
membrane and higher for the system with the an- 
ion-selective membrane [9,10]. Although the mem- 
brane potentials were relatively small due to the 
presence of a background electrolyte [10], i.e., 30 
mmol/L of KCl, the redistribution of [K]t and [C1]t 
in accordance with the membrane potential led to 
the significant change in [Ca]T/[P]T as mentioned 
above. 

For a given cariogenic challenge, i.e., a given 
"plaque" pH, the rate of demineralization is highly 
dependent on the membrane permselectivity 
(Table 4). It is worth noting that although the 
cationic permselective membrane reduces the dif- 
fusion of all anions, it actually led to significantly 
increased fluxes of Ca and P when compared to 
those obtained with the nonselective membrane. 
On the other hand, the anionic permselective 
membrane, which reduces the diffusion of all 
cations (except H"^), decreased the fluxes by a hun- 
dred fold. 

The above findings, which are in good agreement 
with the experimental data, are a direct conse- 
quence of the diffusion-dissolution interactions. 
Previously reported mathematical models for caries 
[3-7,31,32], although possibly more sophisticated in 
some respects, have not been able to demonstrate 
this important phenomenon that has been con- 
firmed in both bench-scale [9-11] and microanalyt- 
ical [1,2] caries models. One model [31,32] utilizes 
coupled diffusion as a basis for describing the ef- 
fects of diffusion on dissolution. This model, which 



becomes mathematically complicated when the 
number of components is three or higher, was not 
developed to study the effects of membrane 
permselectivity on rate of demineralization. 

In this initial assessment of the simulation 
model, the same "plaque" solutions as those em- 
ployed in the bench-scale model were used. In fu- 
ture studies, "plaque" solutions with compositions 
that mimic more closely those of real plaque under 
cariogenic conditions should be used. The effect of 
weak acids such as lactic or acetic acid in the 
"plaque" compartment will differ from that of HCl 
because the weak acid will be only partially dissoci- 
ated in the pH range associated with cariogenic 
conditions. The presence of fluoride in the 
"plaque" solution or the presence of fluorapatite 
(FAp) in the enamel "lesion" should produce sig- 
nificant effects on the composition of the lesion 
and the rate of demineralization [33]. 

Since the demineralization occurs primarily at 
the advancing front of the lesion, the demineraliz- 
ing acid must be transported from the plaque 
through the diffusion barrier (the intact layer and 
the body of the lesion) to the advancing front in 
order to effect demineralization. Similarly, the sol- 
ubilized ions, e.g., Ca^*, H2PO4, etc., must be 
transported from the advancing front of the "le- 
sion" towards the plaque. Thus, concentration gra- 
dients of all the diffusion species exist within the 
barrier. Although the present model is an approxi- 
mation of the diffusion process taking place during 
subsurface demineralization, it clearly demon- 
strated the results of the interactions between the 
diffusion and dissolution processes. A model that 
allows for a diffusion barrier with a finite thickness 
such that concentration gradients of the diffusing 
species may exist within the barrier should further 
improve the ability to predict the effects of various 
factors on the caries process. 

5. Appendix A 

Figure 2 shows the solubility isotherms of OHAp 
expressed as the logarithms of the Ca and P 
concentrations of a series of saturated solutions. 
The curve labeled [0] (Insert) represents the 
isotherm in the 3-component system, Ca(0H)2- 
H3PO4-H2O, where no "foreign" (noncalcium, non- 
phosphate) ions are present. U{-) [34] represents 
the concentration of a foreign acidic component 
such as Cr. Similarly U{ + ) denotes the concen- 
tration of Na* if NaOH is a foreign component. 
Thus the label "-[3]" (" + [3]") is to be inter- 
preted as the sum of the dissociated foreign acid 
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(base) ions and charged ion-pairs at a concentra- 
tion of 10"^ eq/L, i.e., pt/( — ) = 3. If several foreign 
components t/(±) = |XiZ,C,| represents a net 
"fourth component," where C, is the concentration 
of the jth "foreign" ion or ion-pair, z, is its valence, 
and pU{±) =— logf/(±). U{ + ) denotes an ex- 
cess of foreign base over acid in the summation and 
contrariwise for U( — ). Each curve in Fig. 2 repre- 
sents the applicable solubility of OHAp for a given 
pU value. 
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